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Abstract. We consider the Euler equations for compressible fluids in a nozzle 
whose cross-section is variable and may contain discontinuities. We view these 
equations as a hyperbolic system in nonconservative form and investigate weak 
solutions in the sense of Dal Maso, LeFloch, and Murat. Observing that the 
entropy equality has a fully conservative form, we derive a minimum entropy 
principle satisfied by entropy solutions. We then establish the stability of a 
class of numerical approximations for this system. 



1. Introduction 

Compressible flows in a nozzle with variable cross-section a = a(x) are described 
by the Euler equations, which under symmetry assumptions take the following form 

m 

dt(ap) + d x (apu) = 0, 

dt(apu) + d x (apu 2 ) + ad x p = 0, (1.1) 
dt(ape) + d x (au(pe + p)) = 0, x G E, t > 0, 
where the main unknowns are the fluid velocity u and the thermodynamic variables 
e and p. Here, e = e+u 2 jl is the total energy, and we also define the specific volume 
v = l/p and the pressure p, the temperature T, and the specific entropy S. These 
variables are related via the equation of state of the fluid under consideration, for 
instance in the form 

p=p(p, S). 

Provided p is a monotone increasing function of p (for fixed S), the system above 
is a hyperbolic system of balance law with variable coefficients. 



To deal with discontinuous cross section, following [19 we supplement (1.1 1 with 
the "trivial" equation 

d t a = 0, (1.2) 
so that the whole set of equations can be written as a hyperbolic system in non- 
conservative form 

d t U + A{U) d x U = 0. (1.3) 
In consequence, at least within the regime where the system is strictly hyperbolic, 
the theory of such systems developed by Dal Maso, LeFloch, and Murat |7] (and 
also [2TJ HU [20]) applies, and provide the existence of entropy solutions to the 
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Riemann problem (a single discontinuity separating two constant states as an initial 
data), as well as to the Cauchy problem (for solution with sufficiently small total 
variation) . More recently, LeFloch and Thanh [2U HSj solved the Riemann problem 
for arbitrary data, including the regime where the system fails to be globally strict 
hyperbolicity (i.e., the resonant case). The Riemann problem was also solved by a 
different approach by Andrianov and Warnecke pQ. For earlier work on resonant 
systems, see also [T5j |T2l |5] . 

In the present paper, we pursue the analysis of the Eulcr equations in a general 
nozzle, and establish several properties of solutions. We derive the entropy inequal- 
ity that must be satisfied by weak solutions. The entropy inequality is found to 
have a fully conservative form, so that the notion in nonconservative product is not 
needed to state the entropy inequality. In turn, we can obtain a generalization of 
the so-called minimum entropy principle, originally established by Tadmor |26] (see 
also [25]) for plane symmetric fluid flows. In particular, it follows that the specific 
entropy is non-increasing in time. 

These properties are important as far as the stability of numerical schemes is 
concerned, and one of our main results is a proof that a scheme proposed by Kroner 
and Thanh [15 satisfies a variant of the minimum entropy principle. Recall that the 
discretization of systems of balance laws is particularly delicate, and was addressed 
by many authors jTUJ [3J HI |SJ GO [2] . We show that the scheme under consideration 
not only preserve equilibrium states, but also preserves the positivity of the density 
and the minimum entropy principle. 

2. Entropy inequality for nozzle flows 
Consider the system ( |1.1[ ) supplemented with a given equation of state p 



p(p,e) 7 and set U = {p.pu,pe). The flux function and the right-hand side in (1.1) 
can be expressed as functions of U: 

(pu, pu 2 + p(p, e),u(pe + p(p, e))) := f(U), 

(0,p(p,e),0) := g(U). 

The system under consideration can therefore be written as 

d t (aU) + d x (af(U))=g(U)^. 

Recall that weak solutions to this nonconservative systems are defined in the 
sense of Dal Maso, LeFloch, and Murat |7]. On the other hand, as we will see 
later on, the (mathematical) entropy inequality associated with this system has a 
conservative form and does make sense in the framework of distributions. 

Consider first the one-dimensional gas dynamics equations corresponding to a 
constant function a (x <G K, t > 0) 

d t p+d x (pu) = 0, 

d t (pu) + d x (pu 2 + p) = 0, (2.1) 

d t (pe) + d x (u(pe + p)) = 0. 

On one hand, the notion of entropy is motivated from physics, and the physical 
entropy is U — pS, where S is the specific entropy. On the other hand, as was shown 
by Harten et al [11] . necessary and sufficient conditions for a twice differentiable 
function IA C of the form 

U c = pg(S), (2.2) 
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to be an entropy of the usual gas dynamics equations is that g(S) satisfies the 
following properties: 

(i) g(S) is strictly decreasing as function of S; 

(ii) g(S) is strictly convex as function of (1/p, e). 

Moreover, the system under consideration is strictly hyperbolic if and only if it 
admits an entropy of the form (2.2 1. 

Consider the Navier-Stokes equations describing a viscous fluid flow in a nozzle 
with smooth area function a v = a v {x) 

d t {a v p v ) + d x (a v p v u v ) = 0, 

d t (a v p v u v ) + d x (a v {p v u 2 v + p v )) = p v d x a v + v d x (b v d x u v ) , (2.3) 
d t (a ) + d x {a v u v {p v e v +p v )) = vd x (b„u v d x u v ), 

where b v = b v (x) > is given and v denotes the viscosity coefficient. 

We consider the limit as v tends to zero. For simplicity, we drop the subscript v 
and derive the equation for the specific entropy. To this end, we assume that the 
internal energy is given by an equation of state e = e(p, S). On one hand, thanks to 
the equation of conservation of mass of (2.3 1, the equation of momentum in (2.3) 
can be written as 

ap(u t + uu x ) +p x a - v (bu x ) x = 0. (2.4) 
On the other hand, the equation of energy in (2.3 1 can be written as 

ape t + e((ap) t + (apu) x ) + apue x + (aup) x = v (buu x )) x . (2.5) 

The second term on the left-hand side of (2.5) is equal to zero due to the conser- 
vation of mass. Using the thermodynamical identity de = TdS — pdv and v = 1/p, 
we can re- write the equation (2.5 1 as 
dp 

apT(S t + uS x ) H (p t + up x ) + (aup) x + apu(u t + uu x ) = v (b u uu x )) x . 

P 

Or, after arranging terms, we obtain the equation of energy 
P 

apT(S t + uS x ) + - ((ap)t + (aup) x ) + u(ap(u t + uu x ) + ap x - v (bu x ) x ) 

P (2.6) 

= ubui. 



The second and the third term on the left-hand side of ( 2.6 1 are equal to zero by the 
conservation of mass of (2.3) and the momentum equation (2.4 1. Thus, we deduce 

itropy S shoi 

d t S + ud x S 



from (2.6) that the specific entropy S should satisfy 

bv r, 



apT x 

Let g(S) be any smooth function of S satisfying g'(S) < 0. Multiplying the above 
equation by apg'(S), we obtain 

apd t g(S) +apud x g(S) = v apg (S) {~^f u l)- 



Multiplying the conservation of mass of (2.3) by g(S) and then summing up with 
the above equation, we find 

d t (apg(S)) + d x (apug(S)) = vpg'{S)(^u 
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It is clear that if the system (2.3 1 admits a sequence of smooth solution U v 



uniformly bounded in amplitude and converging almost everywhere to a limit U — 
(p, pu, pe) when v tends to zero, then, the function U satisfies the entropy inequality 



d t (apg(S))+d x (apug(S))<0, 



(2.7) 



where g is any function satisfying the above items (i)-(ii). Recall that a weak 



solution of (1.1 1 satisfying the entropy inequality (2.7 1 in the distributional sense 



is called an entropy solution. 



3. Minimum entropy principle for nozzle flows 
We always assume that the fluid is in local thermodynamic equilibrium, so that: 

The function (v, e) > e(v, S), is strictly convex. 

This assumption is equivalent to the requirement that the function (v, e) i— > S(v, e) 
is strictly concave. 

Thanks to the divergence form of the entropy inequality (2.7 1, we can establish a 
minimum entropy principle. We begin with the entropy inequality in a generalized 
form and check that the entropy inequality for (1.1 )-( 1.2 1 coincides with (2.7) for 
the entropy pair 



(pi, J 7 ) = (apg(S),apug(S)), 

where the functions g satisfy the assumptions (i)-(ii) in Section [2] 
Consider the hyperbolic system in nonconservative form 

d t U + A(U)d x U = 0. 



(3.1) 



(3.2) 



The entropy inequality for (3.2 1 has the form 

d t U{U)+ \DuU{U) A{U{.,t))d x U(.,t) 



<o, 



where is a given Lipschitz family of paths, and U is a convex function satisfying 

DfjU(U) A{U) = A(U) T DjjU(U). 

Basic properties of the nonconservative product imply that if there exists a function 
T such that 



D U UA(U) = D U F(U), 



then the nonconservative product 



D V U{U) A(U(.,t))d x U(.,t) 



usual one in divergence form, and is independent of the path 
the entropy inequality takes the divergence form 



d t U{U) + d x T{U) < 



(3.3) 

reduces to the 
Consequently, 

(3.4) 



in the sense of distributions. 



We now check that the entropy inequality for ( 1.1 1 can be reduced to the diver- 
gence form (3.4) for all entropy pairs of the form (3.1 1. This will establish (3.3). 
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The system (1.1 1-( 1.2 ) can be written in the nonconservative form (3.2), with 
U = (ap, apu, ape, a) :— (wi,w 2 , w 3 , w 4 ) Replacing the expression of U in the sys- 
tem (1.1|-(1.2), we obtain 

d t wi + d x w 2 = 0, 

d t w 2 + —d x w 2 + w 2 d x ( — ) + w 4 d x p = 0, 

Wi Wi 

d t w 3 + w 2 d x (—) + — d x w 2 + d x (aup) = 0, 

W\ W\ 

d t w 4 = 0. 

After a tedious but straightforward calculation, we arrive at the following system 
d t wi + d x w 2 = 0, 



d t w 2 + (p p 



p e W 3 W 4 p e W 2 W 4 



2 

w 2 



a , f 2w 2 

OxWl + 

V wi 



PeW 2 W 4 \ 



ll<7 



d x w 2 



PeW 4 „ PpW X 

o x w 3 d x w 4 = 0, 



'"'1 



U'4 



(w 2 w 4 (p p p e w 3 p e w^ p\ W 2 W 3 \ 

w 3 + 2- H 3 2- d x Wx 

V w\ \w 4 wf w? wi I wi / 



( 



w 3 + pw 4 p f w 2 w 4 



w 2 



d x w 2 + (- 

. w 2 ( p p wx \ 

H [P — )OxW 4 = 0, 

Wi V wa / 



w 2 w 4 p + w 2 
w 2 Wi 



d x w 3 



Wi \ w 4 
d t w 4 = 0, x € R, t > 



This system has the canonical form (3.2 1, where the matrix A(U) 
given by: 



an 



0, 012 = 1, ai3 = 0, flu = 0, 



f(tf)) is 



a 2 i = P P 
2w 2 



p e w 3 w 4 p e w^w 4 



a 2 2 

023 

a 2 4 
a 3 i 



P 



PeW 4 _ 
Wi 
PpWl _ 

w 4 

W 2 W 4 / Pp 

Wi \w 4 



H'T 



p e w 2 w 4 



W7 



2u 



H (u -e), 

P 



P 



-PpP> 



IU3+PIU4 p e w 2 w 4 
a 32 = o — = e 



PeW 3 PeW 2 _ p_ 

wf wf w%. 

P _ p e U 2 
P P 



W 2 W 3 



W7 



u{p P 



_ p e W 2 W 4 
°33 — 2 



w 



a 34 



w 2 



Wi Wi 

p + w 2 _ p e u 

Wl p 
PpWl 



u, 



P - 



W 4 



= u(p-p p p), 



p e (u 2 - e) -p 



O41 = CI42 — 043 — 044 — 0. 



These coefficients will be needed in the proof of the forthcoming theorem. 



(3.5) 
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Proposition 3.1. Consider the system \1. ty - §Tjfy in the form (3.Hfy , (3.5). Let g 

be any function satisfying the hypotheses i) and ii) in Section 3. Then the function 
hi = apg(S) of the conservative variables [ap,apu,ape,a) is convex. Moreover, it 
satisfies 

B V U A{U) = D V T{U), T{U) = apug(S), 
which implies that (hi, J 7 ) is an entropy-pair of the Euler system. Consequently, the 
entropy inequality in the sense of nonconservative products can be written in the 
divergence form 

{apg{S)) t + (apug(S)) x < 0. (3.6) 

Proof. First, as is shown in [TT], the function pg(S) is convex in the variable 
(p,pu,pe). Therefore, the function apg(S) — pg(S), p :— ap, is convex in the 
variable (p, pu, pe). Since this function pg(S) can be seen as dependent only on the 
first three variables (p, pu, pe), it can thus be seen as independent of a. Therefore, 
it is also convex in the variable (p, pu, pe, a) = (ap, apu, ape, a). 

Second, the equation of state for the specific entropy being written as S — S(p, e), 
a straightforward calculation shows that 

U(U)=apg{S)=w 1 g(s( W ± ,^ -~(^ f)), F{V)=vU{U). (3.7) 

Using the thermodynamic identity de — TdS — pdv = TdS + jjidp, we have 

S -- S - ^ 

Note also that 



Therefore, it follows that 

■ S p , S e ,w 2 

S W4 

1 ,,„s „ W 



D V U(U) = (g(S) + + -3<?T ~ W3 ))' 

V x w^ w{ 2 

~g'{S)S t l -,g l {S)S t ,-g l {S)S p ^\) (3.8) 



= (9(S) + —^{--- e + u ), — — ) 

and 

D V T{U) = uDuU{U) +U(U)Du( — ) 

^Wl ' 

= (uhi Wl - —^U(U),uU W2 + ,vU W3 ,vLi W4 ). 
wf Wi 



(3.9) 



From (3.5), (3.8 1, and (3.9), we claim that 

B := D V U A(U) - D V T{U) = 0. (3.10) 
Actually, setting B = (61,62,63,64), we have 

61 = -~wig'{S)S t — a 2 i + wig'(S)—a 3 i - uU Wl + —^U 

Wl Wl w{ 

= g'(S)S £ (a 3 i - — a 2 i) - uU Wl + ^hi 
wi wf 

= g {S)S e { 2 1 3 2~) ~ uMwi H 2"' 

wi w? wi wi 
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thus 



61 = g l (S)SJu 3 — ue) — uU w , H %U 

P w{ 

_9'(S)u 2 p g'(S) -p _> 



T 

0, 



V - ~ e) - u(g(S) + + u 2 - e)) + ug(S) 



and 



T p T p T p p 

- {uU W2 H ) 

Wx 

g'(S)u 2 U 
= 9{S) Tf, vUwi =0, 

and finally 

g'(S)up e g'(S) p e u 

h = jf, + +u)- vU W3 

T p T p 

g'(S)u g'(S) 
04 = Ppp—j, h "CP _ PpP) —f 'MAw* = °- 

From the above relations we easily check (3.10), which completes the proof of the 
proposition. □ 

We are now in a position to establish the minimum entropy principle for gas 
flows in a nozzle. 



Theorem 3.2. If U is a bounded entropy solution to the system (1.3), then it 
satisfies the minimum entropy principle." 



inf S(x,t) > inf S(x,0). 

\X\<R \x\<R+t\\u\\ L ao 



We will need: 



Lemma 3.3. Given a real p > 1, consider the function g(S) := (So — S) p , where 
Sq is a constant such that Sq — S > for all S in the domain under consideration. 
Then: 

(i) g(S) is strictly decreasing and strictly convex as a function of S, 

(ii) g(S) is strictly convex as a function of(v,e). 

Proof. We have g'(S) = -p(S - S)?- 1 < and g"(S) = p(p - 1)(S - Sy- 2 > 0, 
so (i) follow imediately. 

Next, since the function S(v, e) is strictly concave as a function of (v, e) for 
< A < 1, for («x, ei) (v 2 , £2) we have 

S(\(vx,ei) + (1 - A)(«a, e 2 )) > XS( Vl , ei) + (1 - X)S(v 2 ,e 2 ). 

Thus, by (i) it follows 

g(S(\(vi, ei ) + (1 - A)(«a, e 2 ))) < g(XS{v lt e x ) + (1 - A)5(« 2 , c a )) 

<Afl(5(«i,e 1 )) + (l-A)s(5(«2,e2)) 1 

which establishes (ii). □ 
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Proof of Theorem \3.S\ Let g — g(S), where S is the specific entropy, be any func- 
tion satisfying the conditions (i) and (ii) stated earlier. We claim that any bounded 
entropy solution of the system ( |1.1| satisfies 

p(x,t)g(S(x,t))dx < / p(x,0)g(S(x,0))dx, (3.11) 

\x\<R J\x\<R+t\\u\\ LX , 

for any nonnegative decreasing function g(S). 

We follow the arguments in the proof of Lemma 3.1 of [25] . We integrate 
the entropy inequality in the divergence form (3.6) over the truncated cone C — 
{(x,t)\\x\ < R+ (t — t)\\u\\loo, < r < t}. Denoting by (n x ,n t ) the unit outer 
normal of C, Green's formula yields 

apg(S)(n t + un x )ds < 0. 

dc 

The integrals over the top and bottom lines of dC give the difference between the 
left- and the right-hand sides of (3.11 1. It follows from the last inequality that this 
term is bounded from above by 

apg(S) (n t + un x ) ds. 

mantle of C 

We will show that the last quantity is non-positive. Indeed, on the mantle we have 

K,n f ) = (i + IHII^)- 1/2 (x/N,|H| L ^), 

thus 



n t + un x = (1 + \\u\\ 2 Lao y 1/2 (\\u\\ L ^ + ux/\x\j 



> 0. 



By the condition g(S) > 0, we obtain the desired conclusion. 

Now, consider the family of function g(S) = (S — S) p , where p > 1 and S is a 
constant satisfying S + So > 0; for instance, So = HSIIlcc + 1- To proceed, we use 
Lemma 3.3 The inequality (3.111 yields 



(3.12) 



p(x,t)(S - S(x,t)) p dx< / p{x,0)(S - S(x,0)) p dx, 

\x\<R J\x\<R+t\\v\\ L °° 
or 

( / p(x,t)(S -S(x,t))Pdx) 1/P < ( f P (x,0)(S a ~S(x,0)) p dx) 1/P . 

K J\x\<B> ' K J\x\<R+t\\u\\ L °a ' 

This means that 

||//"(. > *)(5o-S(. > *))|Up [ _fl^ 

< ||p 1 ^(.,0)(S -S(.,0))|| iPhfl _ t || u | Uoo , i?+t || u | Uoo] . 
Letting p — > +oo in the above inequality we obtain 

\\Sq - S\\ L ^[^ RtR ] < \\S - S , || L =o[_ J? _ t || u || LOc ij . +t || tl || l/00 ]. 
By definition, we can write this in the form 

sup (Sq — S(x, t)) < sup (So — S(x, 0)), 

|x|<_R |a;|< l^oo 

or, suppressing the large constant So, 

S + sup (-S(x,t)) < S + sup (-S(x,Q)). 

\x\<R \x\<R+t\\u\\ L ^ 



THE MINIMUM ENTROPY PRINCIPLE FOR COMPRESSIBLE FLUID FLOWS 



9 



Eliminating So and using (for instance) sup (—S(x,t)) = —mi S(x,t), we arrive at 
the desired result and the proof of Theorem |3.2| is completed. □ 



4. An entropy STABLE AND WELL-BALANCED SCHEME FOR FLUID FLOWS IN A 

NOZZLE 

4.1. Equilibrium states and admissibility criterion. In this section we inves- 
tigate various properties of approximate solutions generated by a finite difference 
scheme for the Euler equations in a nozzle; this scheme was first proposed in [141,115) . 
For definitcness and clarity in the presentation, we consider stiffened gases described 
by 

p= (7- l)p(e-eoo) -7Poo, 1 < 7 < 5/3, (4.1) 
where eac,Pac ar e constants depending on the material under consideration with 
Poo > 0. 

One key property of the well-balanced scheme under consideration is that it 
preserves equilibrium states. In this subsection, we will recall some basic facts and 
explain our selection criterion in the construction of the right-hand state that can 
be connected to a given left-hand state by a stationary wave. 

Let us observe first that the system ( |1.1[ )-( [0| is non-strictly hyperbolic. More 
precisely, the phase domain is divided into three sub-domains so that in each of 
these domains the system is strictly hyperbolic, and along the phase boundary the 
characteristic fields coincide. We express here all other thermodynamics variables 
in term of (p, S). For smooth solutions, the system ( | 1 . 1 [ )- ( p~72] > is equivalent to 

pu 



p t + up x + pu x H a x = 0, 

a 

u t + ^-p x + uu x + —S x = 

P P 
S t + uS x = 0, 

a t = Q. 



(4.2) 



Thus, in the variable U — (p, u, S, a) the system (1.1 1-( 1.2 1 for smooth flows can be 
written in the nonconservative form Ut + A(U)U X — 0, where 



A(U) 



( u 

1± 

P 


Vo 





Ps_ 

P 
u 





up\ 
a 



/ 



The matrix A(U) admits four real eigenvalues, provided p p (p, S) > 0. Therefore, 
for p p (p, S) > the system ( 1.1 1-( 1.2 ) has four characteristic fields associated with 
the eigenvalues 



A = 0, Ai = u - \Jp p (p, S), A 2 = u, A 3 
The phase space is decomposed in several regions: 



y/p P (p, s )- 





= {u 


G 2 


= {u 


G.3 


= {U 


G4 


= {u 



Xi(U) < \ 2 (U) < X 3 (U) < X (U)}, 
X 1 (U) < X 2 (U) < X (U) < X 3 (U)}, 
Xi(U) < X {U) < X 2 {U) < X 3 (U)}, 
X (U) < X\(U) < X 2 (U) < X 3 (U)}, 



(4.3) 
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together with isolated surfaces along which the system fails to be strictly hyperbolic: 

£+ = {U : X 1 (U) = X (U)}, 

S = {U : X 2 (U) = \ {U)}, (4.4) 
£_ = {U : X 3 (U) = X (U)}. 

Note that the matrix of the hyperbolic system is diagonalizable on Eo too. (This 
property is relevant for a flow at rest for which u = even on a discontinuity of a.) 

Next, we consider some properties of equilibrium states. The entropy is constant 
across any stationary wave. So, we may talk about states by ignoring the component 
S. Suppose that a left-hand state U— — (p_,u_,a_) is given, where a = a_ is the 
value of the cross-section. A state U+ = (p+,w+,a + ) (with the cross-section a+) 
which can be connected with C7_ via a stationary wave is determined by the system 

S — S — — S~i~^ 

P = p(p, S-), 

[apu] = 0, ( 4 - 5 ) 



= 0. 



Here, h = e + pv is the specific enthalpy which satisfies 



— h(p,S„) = v— P (p,S_). (4.6) 



To solve the system of equations (4.5 ) for p + = p. we find the roots of the equation 



$(p):=(u 2 _+2h(p-,S-))p 2 -2p 2 h(p,S-) = (^^) 2 . (4.7) 

A basic calculation implies that the equation ( |4.7| has a root if and only if 

a+ > a min {U.) := F (4.8) 

V*(Pmax) 

where p max is the (unique) value such that d$ ^ max ^ = 0. Moreover, in this case, 



(4.7 1 has two roots ipi(U—,a+) < (^2(^-7^+), which coincide if and only if a+ = 
a min (U-). 

Next, we address the question of how to select the right states between <px(U_, a + ) 



and ip 2 (U-,a+). First, given a state U- we observe that the last equation of (4.5) 
also determines a stationary curve u — u(p) in the plane (p,u). Hence, the third 
equation of ( |4.5[ ) implies that the component a can be expressed as a function 
a = a(p) of the variable p along this curve. We then postulate the following admis- 
sibility criterion (see [T21 l^] V 

Monotonicity Criterion. Along the stationary curve in the {p,u) -plane and 
between the left- and right-hand states of any stationary wave, the component a 



determined by (4.5 1 and expressed as a function of p must be monotone in p. 



It is checked in \22\ that: 



Lemma 4.1. The Monotonicity Criterion is equivalent to the condition that any 
stationary wave remains in (the closure of) a single phase. 
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4.2. Definition of the well-balanced scheme. The numerical scheme consid- 
ered in the present paper is defined as follows. The mesh-size is chosen to be 
uniform, i.e., xj+i — Xj = Ax = I, and we introduce standard notation 



i > max{| U J| + fipjj^sf)}, 



A *,« 
At = XAx, 

U := (p, pu, pe), f(U) := (pu, (pu 2 + p),u(pe + p)), (4.9) 

= U? - \(g hF (U?, U? +h _) g hF (U?_ h+ , U?)), 

= U?-i,+ + U? + i,- A 
2 2 

where g hF (U, V) is the Lax-Friedrichs numerical flux: 

9 hF (U, V) = -(/(£/) + f(V)) -~(V- U). (4.10) 
The description of the states 

U? +l! _ = (p,pu,pe)] +1> _, U?_ li+ = (p,pu,pe)2_ lt+ (4.11) 
will be given shortly. 

Remark 4.2. Although the mesh-size is chosen to be uniform, the arguments below 
still hold for non-uniform meshes. For example, one may take At = Ainfi e z{|^;+i — 
Xi\}, provided the infimum is not zero. 



In the scheme (4.9), the states 

Uj+i,- = (P'P u >Pe)" +li _, U?_ 1<+ = (p, pu, pe)£_ lf+ 

are defined as follows. 

First, we observe that the entropy is constant across each stationary jump, and 
we set e? +1> _ = e(/>? +lj _, SJ +1 ), e^ +1< _ = e^ +h _ + (u] +1 _) 2 /2, and so on. Then, 
we determine _ from the equations 

n n n n _ n n n 

«,,) 2 « +1 ) 2 (4-12) 

where h denotes the specific enthalpy defined by dh — TdS + vdp, and is expressed 
as a function h = h(p, S) of the density and the specific entropy. In fact, we have 

, , m V P {P,S) 

hp{p,b) = , 

for any fixed S, and the stationary jump always remains in a given region Gi. 
Similarly, we determine + , + from the equations 

n n n n u n — n n n n n n 

« l) 2 « 1 +) 2 ( 4 ' 13 ) 

by requiring that the stationary jump always remains in a given region G{. This 



definition selects the state in the same phase. Observe from (4.13) that we have 
two states, one belongs to the same phase (or region) and one belongs to another 
phase. 
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4.3. Preservation of steady states. Our well-balanced scheme is stable and 
maintains the equilibrium (steady) states |15j . Considering a stationary wave, then 
one has 



(7™ ri n 7/™ — rt™/1™7/" 

a j+iPj+i u j+i - a j Pj U J ' 



«+i) 



Mp"+i) 



It follows from (4.12) and (4.141 that 



n n 

Pj-1,+ 



Py 
Py 



K) 2 



^•+1,- 



i.e. 



and, thus 



U 



3+1 - 



17? , K 



n+l 



3-1,+ 



u n 



J ' 



This means that the scheme generates the true steady states. 

On the other hand, in the special case that a is a constant, then 



Tjn _ jjn jrn 

u j+l.- ~ u i+l> u i-i. 



u 



i-i' 



and our scheme reduces to a rather standard scheme. 

Let us next illustrate via some numerical experiments the property that the 



scheme (4.9 1 preserves steady states. Consider the Ricmann problem for (1.1) 



where we take 7 = 1.4. and the Riemann data are such that the left-hand side state 
is 

Ul = (pL,u L ,p L ,a L ) = (2,0.5, 2 7 ,1), 
and the right-hand side state 

Ur = (p R ,u R ,p R ,a R ) = (2.080717229626240, 0.320402338758170, p Rl 1.5), 

which is chosen so that these two states are connected by a stationary wave: 



U(x,t) 



Ur, 



if x < 
if x > 0. 



(4.14) 



We compare here the approximations of the stationary wave (4.14) given by our 
scheme ( |4.9[ ) and a standard scheme based on a direct discretization of the right- 
hand side (using central difference, forward difference, or backward difference). 
With the latter, stationary waves are not preserved. For dcfinitcncss, we choose 
the Lax-Friedrichs (LF) scheme. 

First, we note that the scheme ( |4.9| with 1000 mesh points gives the approximate 
density, velocity, and pressure which are virtually the ones of the exact stationary 
wave, respectively; see Figure [T] Figure [2] and Figure [3] 

Second, with the Lax-Friedrichs scheme, the density, velocity, and pressure cor- 
responding to the stationary wave under consideration are not well approximated; 
see Figure [1] Figure [5] Figure [6] The solution exhibits large spikes near the dis- 
continuity since the equilibrium states are forced out of their equilibrium positions 
and generate new waves. 
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2.08 
2.06 
2.04 
2.02 
2 



1 .98 1 — 
-0.5 



"exact solution" - 
"1000 point" 



0.5 



Figure 1. Density of the stationary wave (4.14) - proposed 
scheme with 1000 mesh points 



0.5 



0.45 



0.4 



0.35 



0.3 



- "exact solution" 
■ "1000 point" 



-0.5 



0.5 



Figure 2. Velocity of the stationary wave (4.14) - proposed 
scheme with 1000 mesh points 



4.4. Positivity and minimum entropy principle. Consider stiffened gases, 
where the local sound speed c := \Jp p {p, S) is real. We show that if the initial den- 
sity if non-negative, then the density of the approximate solution is non-negative. 
We also establish the minimum entropy principle for approximate solutions. 
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2.8 



2.75 



2.7 



2.65 



"exact solution" 
"1000 point" 



-0.5 



0.5 



Figure 3. Pressure of the stationary wave (4. 14 1 - proposed 
scheme with 1000 mesh points 
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2.4 
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2.1 
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1.7 

-0 
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"2000 point" 
"4000 point" 



0.5 



Figure 4. Density of the stationary wave (4. 14 1 - LF scheme with 
500, 1000, 2000, and 4000 mesh points 



Theorem 4.3. The scheme (4-9) satisfies the following : 

1. Positivity property: If pj > for all j , then pj > for all j and n. 

2. Minimum entropy principle: for all j and n 

> min^^i}- 
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Figure 5. Velocity of the stationary wave (4.14 1 - LF scheme 
with 500, 1000, 2000, and 4000 mesh points 




exact solutino" 
500 point" 
1000 point" 
2000 point" 
4000 point" 



Figure 6. Pressure of the stationary wave ( 4.14 1 - LF scheme 
with 500, 1000, 2000, and 4000 mesh points 



Proof. 1. We need only show that for any given integer n, if > for j, then 



> for all j. Indeed, the scheme (|4.9| provides the density 



PT = 



> 



n n 4- n 

Pj-1,+ + P. 



i+i, 



> 



p7- 


1,4 


2 
+ 


Pj+l- 


p1- 


1,4 


2 

+ 


Pj+1- 



- - mBx{\^_ lt+ \,\y!j+i,-\W-i,+ + P?+i,-) (4 ' 15) 



(l-Amax(K_ 1)+ |,|< +1) _|)). 
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Since stationary waves provided by (4.12) and (4.141 always connect states with 
dties, it follows from the inequality (4.151 that 

>0 whenever 1 - Amaxflti?]. J, _|} > 0. (4.16) 



Besides, the hypothesis implies that the function h in (4.131 is concave in p: 



h pp (p,S) = ( 7 -2) 



P P (p, S) 
P 2 



< 0. 



It follows from (4.12) and the concavity of the function h that 
1 



;K+i,-) 2 



i 



<-h p {p? +1 ,S? +1 )(p? +li _-$ +l ) 



Pj+1 

< p P {Pj+i,Sj +1 ). 



Therefore, we get 



2p p (pJ +1 ,Sf +1 )< 



Similarly, 



(4.17) 
(4.18) 



From <\<L9j, ( |4.16| , ( |4.17| , and ( |4.18[ ), we conclude that 

p] +1 > for alii, n, 
which establishes the first statement of the theorem. 

2. The proof of the second statement is based on the following classical result: 
assume that U is a strictly convex function in R N , and that there exists a function 
T and a vector-valued map / such that D3- = DU Df. If U is a vector defined by 

V + W A 



U : 



r(/00-/(W)), 



then 



(4.19) 
(4.20) 



Now, comparing (4.9) and (4.191, we deduce from (4.201 that the scheme (4.9) 



satisfies the numerical entropy inequality 



U{Uf_ h+ )+U{Uf +1; _) x 



+ ^(^-i. + )-^JU-)), 



for any entropy pair of the form (4.1 1. Thus, we have 

V Pi ) <2^ a iPi-i,+3(5j-i) + ajPj+i,-5( 5 j+i)J 

1 ^ / n n n / on \ n n n / nn \ \ 

+ ^{a :j p :j _ ly+ u J _ ly+ g{b j _ 1 ) - dj p j+1 _u j+1 _g{b j+1 )), 



or 

,«+i 
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In view of Theorem 



4.3 



p" is non-negative. Taking g(S) = (So — S) p ,p > 1, 



where So is some constant such that So — S > 0, and recalling Lemma |3.3| we obtain 



i/p 



i/p 



+ ^ +1) -(l-A^ +li _)(5 -5; +1 ) 3 '. 

Thus, we get 

+ ^"+i,-(l-A^ +1 ,_)(^o-^i) p ) 

< + + ^"+i,-d - K+i,-)) 

x maxjS'o - S£_i, So ~ S? +1 }. 
Letting p — > +oo in this inequality we obtain 

So — S"^ 1 — niax{So — S^_i, So — "S^+i} = — min{<5'™_ 1 , S'jYi}- 

or 

which completes the proof. □ 
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